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ABSTRACT 



This paper considers an alternating renewal process 
with time to failure being distributed exponential (X) , 
and with a general distribution for time to repair. 

Bounds and the precision of the bounds for the avail- 
ability function were obtained. 

A computer program was written to solve for the 
availability function and some other quantities when the 
down distribution was a gamma a integer. The 

availability function with down distribution being a mix- 
ture of gamma- functions is considered. 
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I. DESCRIPTION OF ALTERNATING RENEWAL PROCESS 



A. SYSTEM 

The system considered consists of a single unit or 
component which is repaired upon failure and then re- 
turned to operation. 

In the model, the following assumptions are made: 

- Time to failure is a random variable (R.V.) X with 
distribution F (specified below) . 

- Time to perform a repair is a R.V. Y with dis- 
tribution G (specified below) . 

- Repair commences immediately upon failure. 

- Unit is returned to operating state immediately 
upon completion of repair. 

- Unit is as_ good as new after repair. 

- Consecutive operating times between failures (as 
well as repair times) are independent. 

Thus the system can be modeled as an alternating re- 
newal process . 

B. DISTRIBUTIONS 

The time to failure - termed uptime - is throughout 
assumed to be exponentially distributed, with probability 
density function (p.d.f.): 

w L- xt . t > o 

f(t) = < 0 , t < 0 

A > 0 . 
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Three different distributions for time to repair - 
termed downtime - were treated: 

- Gamma distribution. 

- Mixture of two exponential distributions. 

- Mixture of gamma and exponential distribution. 
The respective p.d.f.’s are: 



( 1 . 2 ) 



g(t) = 



B 



a-1 -Bt 



( 1 - 3 ) 



g(t) = 



( 1 . 4 ) 



g(t) = 



(a-1) ! t ” ' e 
0 

' B > 0 

a = 1 , 2 , 3 , • • • 

mf^e ^ + (l-m)S 2 e ^ 2t 

0 

' 3 x , 3 2 > 0 
0 < m < 1. 

(l-m)B 1 e‘ Blt + mBjte' 6 * 1 " 
0 

> 0 

0 < m < 1 . 



t > 0 
t < 0 



t > ° 
t < 0 



t > 0 
t < 0 



C. LAPLACE TRANSFORMS; GENERAL 



Denoting the operating state by 1 and the failed state 
by 0, the objective is to find: 

- P 1 1 Ct ) : The ’’availability” of the system, i.e., 

the probability that the system will be ”up” or operable 
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when called upon at a random time t, given that it starts 
in state 1 at t = 0 . 



The following are also of interest: 

- Mjjft): the expected number of visits to state 1 

in period (0,t), given that the system is in state 1 at 
t = 0. 

- P Ql (t): probability system is in state 1 at time t, 

given it starts in state 0 at t = 0 . 

- E m (T): expected time system is "on" in (0,T). 

It is shown in [1] that 



Taking Laplace transform on both sides of (1.5) and 
(1.6), denoting the transform by circumflex ("hat"), i.e., 



(1.5) 




o 



( 1 . 6 ) 



t 




o 




obtains : 



(1.7) 



M lx (s) = 



f (s) *g(s) 



s(l-f(s)-gts)) 



In similar manner: 



( 1 . 8 ) 



Puts) =—i^Q 



s(l-f (s) .g(s)) 
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( 1 . 9 ) 



_ (l-f(s)) gfs) 

s(l-l(s)«g(s)) 

( 1 . 10 ) 

E(s) = - - '-M . 

s 2 (l-f(s)g(s)) 

Only for some of the distributions considered can the 
inversion be accomplished in a closed form. For the gen- 
eral case with gamma down distribution, a computer program 
was written to give numerical solutions. 

If a digital computer is not available, good approxi- 
mations are of interest. These - together with bounds - 
have been developed for the availability function. 

Some solutions in closed form in the simpler cases 
are provided. 
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II. EXACT SOLUTION - LAPLACE TRANSFORM 



Laplace transforms of the 

density functions considered: 

( 2 . 1 ) 

* \ 

£ C S ) = XTF 

( 2 . 2 ) 




(2.3) 



up -and downtime probability 



Ref. (1.1) 
Ref. (1.2) 



g(s) 

(2.4) 

g(s) 



®i 

e i +s 



m 



pi 

p i +s 




b 2 

B 2 + s 




Ref. (1.3) 



Ref. (1.4) 



In the case of (2.2) one obtains upon substituting into 
(1.7) to (1.10) and rearranging: 

(2.S) 



Mil ( s ) 



where , 

( 2 . 6 ) 

q(s) 



XB a 

s-q(s) 



(s+er 



a 

xl 

i = l 



0 1 " 1 (s+B) al 



(2.7) 

Puts) 

( 2 . 8 ) 



(s+B) a 
s* q(s) 



Poi(s) 



B a 

s *q( s ) 
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(2.9) 



B nCs) 



_(s+6i“ 

s 2 q(s) 



In cases (2.3) and (2.4), only P xl (s) has been con- 
sidered: 



( 2 . 10 ) 



P.iCs) 

( 2 . 11 ) 

Put*) 



where 



A 

B 

C 



(s+6 j) (s+6 2 ) 

; Ref. (2.3) 

s[s 2 +(X+6+3)s+6(X+6)-Xm(6-3)] 



(s+6j) (s+B 2 ) 2 
s [ s 3 +As 2 +Bs+C] 



; Ref. (2.4) 



+ 2 $ 2 + ^ 

6 2 + 26 2 ( X+8 j) + Xgjm 

X6 2 (1-m) + 6 , 6 2 + 2X6 6 m. 

12 12 



9 



III. INVERSION 



A. COMMENTS 

Consider the transform of a rational function P(s)/sq(s), 

/N A 

where the degree of P(s) is less than sq(s) . This is a suf- 
ficient condition for the use of the residue theorem in 
finding the inverse Laplace transform [2], thus: 

(3.1) 

F (t) = E (Residues of e st f(s) at poles of F(s)). 

As is shown in Appendix B only simple poles exist if 
B > X. This simplifies the inversion. 

Equation (2.9) has a double pole at s = 0 . The ap- 
propriate adjustment to the procedure used is considered 
later. 

In the sequel, primarily the inversion of Pj^s) will 
be considered. However, since the method used is applicable 
to all cases where Laplace transform is of type (3.1), 
the discussion in some detail of inversion of Pjjfs) will 
show how the other results were obtained. 



B. INVERSION FORMULA 

Referring back to (2.7) 

(3.2) 

Pll (t) = L- 1 

\ s q( s ) i 

q(s) is a polynomial of degree a. It is shown in Appendix 
A that all roots are complex, except when a is odd, in 
which case one of the roots is real and negative. 
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Factoring sq(s) into linear terms: 

(3.3) 

sq(s) = sCs-rJ (s-r 2 ) • • • (s-r a ) 



where 



r^ = a^ + b ^ i (b l - 0 if a is odd) 

r • i = a . - b.i j = l,2,***,a-l. 

j+i j j j 



Using the residue theorem: 

P 1 j (t) = Z Residues 
a+1 

lim 



a+1 

-I 



(s-r k )(s + 6) a e st 



k=l 5 ‘k 



s^ ri s(s-rj) (s-r 2 ) • • • (s-r k ) • • • (s-r a ) 



(3.4) 



P n(t) = + 



y <v 

fei-fTr 



a r V t 
B) a e K 



6 + “ x ^T a T r k (r k' r )l ) 



k=l 



H = 1 , 2 , • • • , a If k 



Denoting 

(r k +e3 “ = v k + V 



a 

i [ r k^ r k- r P = x k 

k=l 



y k i *■ t k 



then 



(3.5) 

p i.^ = e4x + ^ (v k x k +v k>’k }cos V + (v VV'., ?' 
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(3.5) Cont 1 d . 



sin b k t) 



.V 



x k + ^k 



k = 1 , 3 , 5 , • • • , a-1 . 



Incidentially , for a bounded solution to exist the 
real part of all roots must be negative. 

The numerical solution is contingent on finding the 
roots of sq(s) . This can be done algebraically for a < 4. 

The method used is outlined in Appendix C. In addition, 
for a = 2 the inversion of (3.5) in explicit form is shown in 
Appendix D to illustrate the method. 

For a polynomial of unrestricted degree the computer 
program described in Section V was written. 

If a computer is not available, approximating formula 
Pjjft) are given in Section IV. 

For E 1X (T), equation (2.9), the double root at s = 0 
leads to: 



E „m 



lim — eST 

s->0 ds s 2 (s-r x ) (s-r 2 ) • • • (s-r^) 



a 



i 

lr= 1 



lim 

s+r. 



(s-r*)(s+e) a e sT 

TT s'Cs-r*) 

Sri 



which upon derivation and letting s -* 0 lead to 



( 3 . 6 ) 



E,x(T) = 



Xg(g+?gX (g-1) ) + Xg T 
(6+gX) 2 g+aX 



I 



lim 



(s -r v ) (s + g) a e sT 



T'i s ^ r v s 2 ( s ’ r JCs-Tj) • • • (s-r,) ■ 



k=l 



' (s -r ) 
v or 



A similar procedure as outlined for P xl (t) is used for 
the last term (modified by multiplication by s) . 
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IV. APPROXIMATIONS AND BOUNDS 

A. GENERAL DOWN DISTRIBUTION 
Consider 

P,,(s) = 

s(l-£(s) *g(s) ) 



where uptime X ~ exponential, downtime Y - G, unspecified. 
Then , 



P.iCs) 

(4.1) 

Pll(S> 

By division: 



1 



X 

X+s 




i-x 

s 



(4.2) 



P(s) = 



-x I z& hl + X 2 P--£(.. S I 
s I s 



3 l-g(s) 



Observing that l-g(s)/s 2 is the Laplace transform of the equi 
librium distribution of G (Reference 4) except for a factor 
1/E [Y] , then 

(4.3) 

P„(s) = 1 - RG e (s) + R 2 (G e (s)) 2 - R 3 (G g (s) ) 3 + • • • 



where 



R " Ifr} 



G e ( s ) 



= 1- Sts) 

E LY] • s 2 
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With R < 1, convergence is ensured. This implies that 
P lx (t) > ^ as t « (steady state result). Inverting (4.3) 

Pu(t) = 1 - RG e (t) + R 2 Gg 2 ^ (t) - R 3 G e (3) (t) +••• 

where denotes the ith fold convolution of G with 

itself. 

Letting t «>, 

lim P n (t) = l- R+R 2 -R 3 + **- = g- pf I* - ] [ 7T 

one obtains the steady state result. 

B. BOUNDS AND PRECISION 

It is a reasonable proposition to try to invert the 
first term in (4.3). Inversion of the other terms may be 
difficult unless G 0 is a simple distribution. 

Even though (4.3) is not a series in the time domain, 
by probabilistic considerations the effect of truncation 
may be deduced. 

Noting that R^ (t) > R^ + ^ + ^ (t) , then the ap- 

proximation of P xl (t) will give a lower or upper bound on 
Pjj(t) depending on whether the last term is negative or 
positive , e . g . : * 

(4.4) 

P xl (t) > 1 - R-G 0 (t) 

(4.5) 

P,,(t) < 1 - R-G e (t) ♦ R 2 Gg 2 ^ (t) . 
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The precision with which the bounds are established 
can be obtained. 

Suppose one truncates after one term. Then the 
maximum value of the remainder is less than R 2 . Thus, 

(4.6) 

1 - RG e (t) + R 2 > P.jft) > 1 - RG e (t). 

Having found the lower bound for P^ft), then this value 
is accurate to within R 2 . In a rather typical case with 
R = 1/10 (i.e., steady state = 0.909), the lower bound 
is established to within one per cent. 

If it is feasible to compute the upper bound by 
taking two terms in the expansion, then - by similar 
argument - the precision of the upper bound will be R 3 , 
i.e., to within 1/10 per cent in the case considered. With 
R relatively large, the precision is poor. If feasible, 
more terms must be computed for a desired precision. 



C. GAMMA DOWN -DISTRIBUTION 



Consider case with Y ~ Gamma (a, 3) a = 1,2,3, •••: 



(4.7) 

/s a* X 

R-Gfs) = 

6 





^ (6+s) a -e a 
s 2 (B+s) a 



Noting that s = 0 is root, one may simplify (4.7) as fol- 
lows : 



(B+s) a -B a = B a +B°‘~ 1 s+B 01 ~ 2 CB + s)s + B ct ~ 3 (B+s) 2 s+- • •+f8+s) a ~ 1 s-B a 
s 2 (B+s) “ s 2 (B + s) a 
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/ 



a - 1 



(4.8) 

RG e (s) 



1 1 
s 8 





a 




+ • • • + 




Substituting into (4.3): 




Each term within the bracket contains gamma trans- 
forms, making the inversion simple. 

As an example, consider case with a = 2. Then (4.9) 
becomes : 



(4.10) 



s-P,,(s) 




The resulting inverse for first approximation: 

(4.11) 

P„(t) > 1 - | (1 - e' Bt + 1 - e‘ 6t (Bt+1) ) . 

Using the correspondence between poisson and gamma dis- 
tributions, the value can be obtained from cumulative pois- 
son tables. 
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D. 



MIXTURE OF GAMMA DOWN- DISTRIBUTION 



Consider a mixture of two gamma-distributions with 
parameters (K,3) , (&,y) and weighted (1-m) and m respec- 
tively. The first term in the expansion becomes: 

(4.12) 




The first term is seen to be the linear combination of the 
two terms in the down-distribution. For case K = 2, 
i - 3 , the inversion is: 



(4.13) 

Pn(t) > 1 - H t (t) + ^ H 2 (t)| 

where 

H j (t) = 1 - e‘ Bt + 1 - e' 6t (6t+l) 

H 2 (t) = 1 - e' Yt + 1 - e‘ Yt (yt+1) + 1 - e‘ Yt * 

+ yt + i | . 

In case more terms are required, the inversion 
will be somewhat complicated in that there will appear 
cross-terms between the two distributions. The precision 
of (4.13) will be 

X> ( li^i 2 + si ] • 
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E. WE I BULL -DISTRIBUTION ; SPECIAL CASE 

A lower bound can be obtained when downtime is dis 
tributed Weibull (3,2). The reliability function is: 

(4.14) 

G(t) = e' 6t2 , 6 > 0. 

The expected value is: 

E[T] = YT . 

The equilibrium distribution: 

G e (t) ■ 



/ «' 6X! 
Em 



which becomes upon re-arranging: 




s ~hx 2 (23) 




(4.15) 



G e (t) =VT / dx. 

o 

Thus, G g (t) is ~ normal (0,1/23). The lower bound becomes 

(4.16) 

p„(t) > i - xVT z 

where Z(t//2]3) is normal (0,1). Thus standard normal 
tables can be used for evaluation. The corresponding 
precision is A 2 • tt / 3 - 
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1 . Weibull -Distribution ; Decreasing Failure Rate 



For the special case that the Weibull distribution 
has shape - parameter K = 1/p, p = 1,2,3,** # , the lower 
bound with one term in expansion can be found. The equi- 
librium distribution: 



(4.17) 



G e Ct) 




/ 



1 

e ’Bz p dz 



which upon substitution becomes: 

(4.18) 

1 
P 

B p -Bx p-1 , 

— e x F dx 

(P-1)! 

where 

z = x p . 2 

Equation (4.18) is the gamma-distribution in x(=z p ) with 
parameters (p,B). The lower bound becomes: 




p„(t) > 1 - A r ( t ;p,B) 



with precision: 



El 



As an example, let k = i.e. 

> 1 - x w 

The value of term in brackets 
tables . 



, p = 2 . Thus , 

- e' 6/F je/t + l) . 
can be obtained from poisson 
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V. COMPUTER PROGRAM 



The computer program was written in FORTRAN IV and 
run on the IBM system/360 model 67 computer at the U. S. 
Naval Postgraduate School. The program listing is shown 
in Appendix E. 

The main sequence in the program is : 

(i) Input data: a, 3, A, time increments and number 

of time steps (named respectively M, BETA, UPRATE, TSTEP , 

L in program) . 

(ii) Generate coefficients of polynomial in denom- 
inator of (2.7). 

(iii) Call library subroutine DPOLRT to solve for 
roots, using Newton- Raphs on ' s method. 

(iv) Compute and store coefficients for cosine, 
sine terms in P xl (t), P 01 (t), M xl (t) and E xl (t); e.g., 
see (3.5) and (3.6) 

(v) Sum cosine coefficients for P 11 (t) and add to 
steady state value of P (t) , to give an indication of ac- 
curacy of final result. Theoretically the sum should add 
to 1. Flagging value IR provided by subroutine to in- 
dicate accuracy of root extraction*should read zero. 

(vi) Enter DO-loop for incrementing time and for 
final calculations. 

(vii) Print results: t, P lx (t), P 01 (t), M xl (t), 

E n (t), named: T, PROBUP, DPRBUP, RNWAL , EXPUP respec- 

tively in the program. 

Sample output is shown in Appendix F. 
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APPENDIX A: NATURE OF ROOTS 



Polynomial considered: 

(A. 1) 

P(x) = x a + Xfx 01 " 1 + gx a ’ 2 +'**+ B a l ) ; 6 > X . 

Multiply both sides by x-8, to obtain: 

(A. 2) 

Q(x) = x a+1 - (B-X)x a -XB a = P (x) (x-B) . 

Let a be even . 

Then by Descartes’ Rule of Sign [3], with only one 
variation of sign, there exist at most one real, positive 
root (clearly x = 3) . 

Consider 
(A. 3) 

Q(-x) = -|x|“ +1 - (B-X)x a - XB a . 

Since there are no variations in sign, no negative 
real roots exist. CONCLUSION: For a even, (A.l) can only 

have complex roots. 

Let a be odd . 

By (A. 2), one variation indicates the positive root 
x = 6. Consider 
(A. 4) 

Q(-x) = x a+1 + (8-X) |x a | - XB a . 

One variation imply one real, negative root. CONCLUSION: 
For a odd, (A.l) has one negative, real root, (a-1) complex 
roots . 
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The roots of (A.l) are of the form 



r . = a. 



V 



j = l,2,***a with b 1 = 0 if a is odd. 



Since x = s+£, the roots of s will be 



= r j 



- 6 = 



a j 



b . i 
J 



and 



a 

I rj = - A - aB. 
j-1 



That all a^ < 0 can be verified by Hurwitz Criterion (3) . 
However, since E[X+Y] < 00 and F*G(t) is non-lattice, the 
system will be stable, implying that a^ < 0 for all j [4] 
(* denotes convolution) . 
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APPENDIX B: NON-EXISTENCE OF MULTIPLE ROOTS 



Polynomial considered: 

(B.l) 

P(x) = x“ + X ( X *- 1 + Bx“' 2 +”* + b“' 1 ); B > X. 

For simplicity, let a = 6, and assume that a double and a 
single complex root exist: 



r l , 2 = a 4 bi 
r 3 , 4 = a ± bi 
r 5,6 = c ± di. 

Consider the usual relationship between roots and coeffi- 
cients [3] : 



(B . 2) 



4a + 2c = -X 
(B.3) 

(a 2 +b 2 ) 2 (c 2 +d 2 ) = XB 5 
(B • 4) 



; sum one root at a 
time 



; sum 6 at a time 
(i . e . , product of 
roots) 



2 (a 2 +b 2 ) + (c 2 +d 2 ) +4a 2 +8ac = BX ; sum roots two at a 

time . 



Using Cauchy’s inequality on (B.3): 
(B.5) 



2 (a 2 +b 2 ) + (c 2 +d 2 ) 






(a 2 +b 2 ) 2 (c 2 +d 2 ) 



By (B . 4) , then 




5 
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(B.6) 



2 (a 2 +b 2 ) +c 2 +d 2 = 3X - (4a 2 +8ac) 

= 3X - 4 ( (a+c) 2 -c 2 ) < 3X 

since (a+c) 2 -c 2 > 0. 

Equations (B.5) and (B.6) are in contradiction, since 
3^X3 5 > 3X for all values of 3 > X. Thus double roots are 
not possible. The proof can be extended to any a and any 
multiplicity of complex roots considered. The proof holds 
also for a odd. 

The two inequalities established for a = 6 hold in 
general since all terms of form (a 2 +b 2 ) will appear in 
isolation from all possible cross terms of form (a*c) . 
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APPENDIX C: POLYNOMIAL DEGREE 3 AND 4 



Provided 3 > A, the roots of polynomial of degree 3 
and 4 can easily be found [5]. 
a = 3 



P(y) = y 3 + Ay 2 +• A3y + A3 2 ; y = s + 3. 



Denote : 



a = j (3A3-A 2 ) 



b = 27 (2X 3 -9A 2 3+27A3 2 ) 



A = 



B = 



j + VET7H 



4 27 



- -'Gnr)* 

4 27 j 



Then the roots, solving for s, are: 

r 1 = A+ B- 3- -| (real root) 

r 2,3 = ' 1 ( A+B) ± I ( A ' B ) Y3 _ i-6-| (complex) 

a = 4 

P(x) = x 4 + Ax 3 + A3x 2 + A3 2 x + A3 3 ; x = s+3. 

Denote : 

a = | (X-6g) A8 2 
b = jy- (9B-4A)A 2 B 3 
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y = m cos 0 + X3/3 



R = [4 - + y Y 

D 1 2 = | J X 2 - X6 - y i: (4X 2 g-8Xg 2 -X 3 )/4R| \ 

Then, the four complex roots, solving for s, are: 




The simplicity in form of polynomial 
P(x) = x“ +1 - (6-X)x a - Xg a 

make solving for roots for a > 4 by manual methods easy, 
for example, by Graeffe’s Method [3]. 
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APPENDIX D: AVAILABILITY FUNCTION IN EXPLICIT FORM; g = 2 



In order to illustrate the inversion method used, con- 
sider PjjCt) for a = 2: 

(D.l) 



1 


I 


t 


< (s + B) 2 1 


1 = l" 1 < 


(s + B) 2 


] s((s + B) 2 + X(s+B)+XB) 


[ i 


I sfs-rj) (s-r 2 ) 



where 

r i 2 = j ( -X-26±/4A8-X2 i) ; 6 > X 

= | (-X-26±Ai) ; A = / 4 XB-X 2 



Then , 



Pii(t) 



= lim 
s->0 



s ( s + g ) 2 e st + 
s (s-r,) (s-r 2 ) 



+ lim 
s+r 2 



(s-r 2 ) (s+6) 2 e 
s (s - r j) (s -r 2 ) 



lim 

s_>r i 

st 



fs-r (s + g) 2 e st 
s (s -r j) (s-r 2 ) 



(r 1+ B) 2 t 

+ — ? v e 1 

r j r 2 r 1 ( r 1 ~ T z i 



(r 2 + 3) 2 

r 2 ( r 2 *r 1 ) 



e 



r 2 t 



Since the last two terms are complex conjugate numbers, 
they may be collected and expressed in terms of trigono- 
metric functions: 



p nW 



6 2 

+ e 

B 2+ 2XB 

A t 

(r 1 + B) 2 r 2 (r 2 - r 1 )e 2 + (r 2 +B) 2 r t (r x -r 2 ) e 

r i r 2^ r i‘ r 2^ ^ r 2~ r i^ 
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M> 



$ -(3+y]t (ac+db)cosyt + (ad-bc)sinyt 

+ 2e z z z 



3 + 2 X 



2 +d 2 



where 

(rj+3] 2 = a + bi ; (r 2 +3) 2 = a - bi 

r i( r i- r 2 ] = c + di ; r 2 (r 2 - r x ) = c - di . 

Substituting in actual values, the result is: 

= 3 + tr + A(3 + 2 A) | A cos 2 Z + sin jtj’ e 



-CB+f)t 



Closed form solutions were also found for: 



lll( t) = - + _*jL_ t + A 



Ae 



(64n 



(3 + 2 A) 2 3 + 2A (3 + 2 A) 2 A 



C 2 3 2 + A 2 ) sin jt 



+ A* (23+A) cos yt . 

P ..W = bTII • IFTIiA (C26 + X)sin f + A cos ftj e 



-(e4)t 
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APPENDIX E: COMPUTER LISTING 



C 



1000 



15C0 

C 



2001 

2002 

2003 



I M Pt TCTT o^AL* 8( A-H t O-Y) ,C0MPl?X*16(Z) 

DTMFNS 1 ON Xrv? < 100 ) , CC?F( 100 ) , ROOTR ( 100) , ROOT I ( 10 0) 

DT MENS’ ON Z*M 100) , ZROOT ( 100) , ZNUM(IOO) , ZDNMdOO ) 

DIMENSION SCOE c ( 10C),QC0EF( 100) 

DIMENSION O-FP(lOO) 

D T M f NS T ON DQCOr F { 100) ,DSC0FF(100) ,0QEF C (100) 

DIMENSION ROOPF( 100 ) • P.QC0EF ( 100) ,RSCOEF( IOC) 

01 ME NS! ON ^COO^FdOO) , SSCOE F ( ICO > , EOEF*= ( 100 ) 

DIMENSION XMUM ( 2) , XDN M( 2) » XRDNM ( 2) , X° (2) 

EQUIVALENCE ( Z ZDN'M , X PN M ) , ( ZZNUM, XNUM) , ( ZP ONM, XRDNM) , ( ZZR , XR ) 
M = 5 

U°RATE=1. 0^0/4.300 
P=TA=1C .000 
T STc D = C. 01 P C 

GENERATE POLYNOMIAL AND SOLVE FOR ROOTS 



TU D =0 

J=M+1 

XCOF( J )=1.0DO 
XCOF( V) =UPR ATF 
S'JMl=0.Cr0 
no 1000 M=2.M 
I U?=I Un+1 
L=M-n +1 

XCPF( L )=UPR ATE *BETA**I UP 
CONTI NIUE 

mow cnp INVERSION OF LAPLAGP TRANSFORM 
CALL D- r LRT ( xfco c * C0 C » ^ * POOTR ? p GOT I . TFR ) 

DO 1500 T =1 

ZPOPT ( T ) =HC MR L X ( D OOTP ( I ) •POGTI ( I ) ) 

IQ ( T ) = 7 prrT ( T ) -BET A 
ZNU V ( I ) =ZROOT( I ) 
f GNT T MU P 

' “AB OVF DID CALCULATE NUMERATORS. NOW FOR DENUMERATORS 
DO 2000 T = 1 , m 
tc(T . M ~ .1 ) O n TO 2001 
Z ° R 0 0= Z R ( 1 ) - Z R ( 2 ) 

IF(M.GT.R) r '0 T- p 2002 
ZDNM(1 )=ZP f ? rn -ZR(l) 
on T n ^502 
Z D R 0 0= Z F ( I )-ZR( 1 ) 

1 SR0P=2 

IF(M.E0.2) GO TO. 2501 

GO T o 2003 

0YM‘1Y = C • 0 

IPR0n=3 

ni j M =0 • 0 

DO 25 00 -J- 7 PR OH • M 

IF( ( J.FO.I J.ANO. (J.EO.M) ) GO T® 2501 
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I<=( J.E0.7) GO TO 2500 

Z D PO!>ZPPO n *(ZP< I )-Z p < J ) ) 

2500 CENT? NL ,n 

2501 ZDNM(I)*ZP®eP*ZR(I) 

2000 CONTTMUE 

r NEXT POLLOWS SINE AND COSINE COEFFICIENTS 

!CHFCK=( M /2)*2- M 
T p ( TC MFC* ) 350 1 * 350 2* 3501 

3502 IMJG=1 

GO m 3503 
3501 I M AG=2 

3503 DIJMMY=1.C 

I L TM =M -1 

r>n 350^ T =I MiG, I L I M * 2 
Z7DN^=ZP^m( t ) 

ZZMUM = ZNUM( I ) 
zpdnm=zpnm< T)^D (1) 

ONUM=XNUM< n*XPNM( 1>+XNUM(2)*XDNM<2) 

ODNM=XPNM( 1 )**2+XDM'M 2 ) * * 2 

oroF = ( T ) =2 . o^o-omum/odnm 

S MUM = XN IJM < l ) * XD N M < 2 ) - X \' UM < 2 ) * XP N M < 1 ) 

SChfc ( t j = 2 .GDO-SNUM/ODNM 

DOrpFP( ! ) = 2 .QD 0*PE TA^^M^XDNMI 1 ) / ODNM 

CSrOFF( I ) = 2 .ODO* Q E TA^M^XDNM ( 2) /ODNM 

RDNM=XRDNM( 1 ) *^2+X RDNM ( 2 ) **2 

ROC OF c ( I ) =2.00 Q^’iPR ATE* PET fc**M*X D PMM f 1 ) / PPNM 
P S CPf c ( T i = 2.^O0i'MP D AT = «RPT =£ X~ n N M { 2 ) / D 0 m m 

- QNUM=XMUM ( { ) * X R^-NmTi ) +XNLM ( 2 ) *XR DNM (25 
ESNU M =XNMM( 1)«x^DNM{ 2)-XMUM( 2)*XRDNM( 1 ) 

c scnr-p ( i )=2 . ooo*FSN!jM/ppMM 

FQCCF F ( t ) =2 . 00 0*F ONLJM/RDNM 
SUM 1= SU M l+ n C n E r ( I ) 

3500 CONTINUE 

C NOW FOR SINGLE REAL COE c FICIENTS 

I c ( I m AG.M*=.2) GO TO 3504 
Z ZDNM=ZDN'M 1 ) 

ZZNUM=ZNUM( 1) 

ZRnNM=ZPNM( 1)*Z D < 1) 

QFPP(1 ) =X Ni I M ( 1 ) / X D N M ( 1 ) 

DOFF c ( 1) =P '.*TA**M/XnMM< 1) 

RQFF C ( 1 )=I»F? ATE-PETA^*M/XRDNM( 1 > 

FOEF^ ( 1 ) =XN U" ( 1 ) /X RDM" ( 1 ) 

S<JM l= SUM l+ n E cc ( 1) 

r PRELIMINARIES TO LOOPING POR TIM- 

3504 n!jMMY=0.0 

PMNPl=P c TA/( BE t *+D FLOAT! M ) *UPRATE ) 

ACUP r Y=PMN01 +SUM1 
WRITE (6,1) 
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1 FODMATCl* *10X, 'CHECK c OR ACCURACY OF NUMERICAL SOLUTION*///) 

WP T TE ( 6 , 2 ) ACIJP C Y, tc R 

2 =orma T( . . ,iox , • TNDTCATOR c OP ACCUR AC Y= • t - 16 . 12, / /, 1 1 X, « SUBRO UTI NE 
1 FLAGGING VALU“=* ,!13/////I 

W° TT r ( 6 t 3 ) 

3 FORMAT ( • 11X , 'TIME', 3X, 'AVAILABILITY' , 5X, ' AVAIL /DOWN* ,6X, 'E(RENS 

1WAL ) • t7Xt f F(TI ON) * //) 

PCNST=UPRAT7«r=TA/( B c TA+OFLOAT( M)*U°RATF ) 

RrONST=UFR/'Tr*nei.. r| A T (M>* (BFTA+O .5D0*UPRATE*( QFLOAT( M)- l.OOO) ) / 

1 ( P5ETA+DFL0AT( M ) *UPRATF ) **2 

ECNST=UPP ATE*DFLOAT( M)*( DFLOAT( M )+l.0D0)/( 2 . ODO* ( QETA+UFR ATE^DFLO 
IAT(M) >**2) 

00 4500 L= 1 , 10° 

T=TSTF D *OFCC AT ( L ) 

RMND=0. 000 

T "RM=0 •CO 

0 D MND=0 .000 

OTEP M = 0. 0 DO 

pPMNn=o .COO 

RTF R M =0 • 000 

EPMNOO.COO 

FT c RM=0 *000 

OG 5000 1 = 1 MAG, I L I M ,2 

ZZR=ZP ( T ) 

OPS=0CGS (X R (2 ) *T ) 

SI N=0 ST N( X~ (?) *T) 

^XRnr’=r'rXR(XR( 1 ) -T ) 

RHN0=R w NO+=y MQW- i QOS-OCOEF ( ! ) +SIN*SC0EF( I ) ) 

DRMMP=DR MMQ4-- XPON* ( QOS*DOCOFF { I ) + S T N* DSC PE F < I ) ) 
r r MM n= p p MN|£ + r x PC m* ( Qns*R0C n 5F ( I ) + S! N*R SC 0 - f ( I ) ) 

FRMND = i P MNO+ r XPD^ (COS*FOCOE c ( I )+S IN*ESCCEF( I ) ) 

5000 CONTINUE 

1 F ( I MAG # f'3 . 2 ) on JO 5001 

R c XPQN=DEX° ( ( R00TR(1 )-BFTA )*T ) 

TERM=OFF c ( 1 )*p«=x=>0N 
DTRPM=noFF = u )*R?XPQN 
RTFRM=RQcF=( 1)*R -X? p N 
pJcpM=rOc F ” ( 1 ) *CC X pnM 

5001 PP ORUP = P v.NOl + P M N r ‘+T r °M 
nop c^np-R win i+o o mmq+DTFPm 

p NH At = -F r n N$T + P ON ST^T +RP MISJD+P TFP M 
»=XPUP = rc.\ , ST + RMM<d*T + !r P^N0+cT C RM 
V ,c T T" ( 6,5003) T , C P OR UP ,0 PR R UP , P NWAL , FX PUP 
5003 format ( » * » 1 OX » F6.3,4X* F8 • 6 » SX*F8«6t7X*F9. 6,9X, C 8.6) 

4500 C 0NT7 MU r 
ST°o 
END 
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APPENDIX F: SAMPLE RESULTS 

CHECK c 0? ACCURACY OF NUMERICAL SOLUTION: 

INDICATOR r CR ACCURACY^ 1.000000000000 
SUBROUTINE FLAGGING V*LUE= 0 



TIME 


AVAILABILITY 


AVAIL /DOWN 


F(RENEWAL) 


F(TJ M F ON) 


0 .010 


C .997780 


C . 000000 


0.000000 


0.009989 


0.020 


0.995565 


0.000002 


0.000000 


0.019956 


0.030 


0.993356 


0.000016 


0.000000 


0.029900 


0.040 


0.991151 


C. 000061 


0.000000 


0.039823 


0.050 


0.988951 


0.000172 


0.000000 


0.049723 


0.060 


0. 986756 


C. 000394 


0.000001 


0.059602 


0.070 


0.984567 


0 .000783 


0.000002 


0.069458 


0.080 


0. 982384 


0.001407 


0.000005 


0.079293 


0.090 


0.980207 


0.002335 


0. 000009 


0.089106 


0.100 


0.978038 


0.003645 


0.000015 


0.098897 


0.110 


0. 975 £77 


Q. 005410 


0.000025 


0.108667 


0.120 


0.973725 


0.00"7706 


0. 000040 


0. 118415 


0.130 


0.971 584 


0.010603 


0.000060 


0. 128141 


0.140 


0. 96^455 


C. 014166 


0 • 00 0087 


0.137 8^-7 


0.150 


0.967339 


0.0 18453 


0.000123 


0.147531 


0 . 160 


0. 965238 


0.023512 


0.000170 


0. 157193 


0.170 


0.963154 


0 • C 29 386 


0.000229 


0. 166835 


0.180 


0 .961089 


0.036106 


0.000301 


0. 176457 


0.190 


0.959044 


0. C43692 


0.000390 


0.186057 


0.200 


0.957021 


0.052157 


0.000496 


0. 195638 


0.210 


n. < 5 5 ^ n ? ^ 


0 .061 504 


0.00C622 


0.205196 


0.220 


C. 553050 


0.071727 


0 * 00 C770 


0.214738 


0.230 


0.951106 


0.082810 


0.000942 


C. 224259 


0.24 0 


C. 549192 


0.094731 


0.001139 


C. 233760 


0.250 


C .947309 


0.1 0"746 1 


0.001363 


0.2^3243 


0.260 


0.945^60 


0.120962 


0.001617 


0.252707 


0.27 0 


C. 543645 


0.135194 


0.001901 


0.262152 


0.280 


0.941867 


0. 150110 


0. 002218 


0.271580 


0.290 


0.94012*7 


0.165660 . 


0. 002569 


0.280990 


0.300 


0. 538425 


0.1 81790 


0.002955 


0.290382 


0.310 


0.936*764 


C .198444 


0.003377 


0.299753 


0.320 


C. 535144 


0.215565 


0.003837 


0.309118 


0.330 


0.933566 


0.233093 


0. 004336 


0.318461 


0.34 0 


0.932031 


0 • 2 5 r, 96 9 


0.004873 


0.327789 


0.350 


C. 530 540 


0.269135 


0.005451 


0 .337102 


C.360 


0.9 29G Q 2 


0.287531 


0.006070 


0.346400 


0.370 


C. 527688 


0.306100 


0.006729 


0.355684 


0.380 


C. 526329 


0.3 247 8 "* 


0.007430 


0.364954 


0.390 


0.925015 


0.343536 


0.C08173 


0.374211 


0.400 


C. 523745 


0.362297 


0.008957 


C. 383454 


0.410 


0.922520 


0.381019 


0. 009783 


0.392686 


0.420 


0.921338 


0.399655 | • 


0.010650 


0.401905 


0.430 


C. 5202C1 


418162 


0.011559 


0.411113 


0 .440 


0.919107 


C. 436497 


0 . C 1 2 5 0 9 


0.420309 


0.450 


0. 918C56 


0.454622 


0.013499 


0 • 42°49 5 


0.460 


0. 917C47 


0.472502 


0.014529 


0.438670 


0.470 


0.916030 


0.490103 


0.015599 


0.447836 


0. 480 


0. 915153 


0.597396 


0.016707 


C .456992 


0.490 


0.914267 


C. 524355 


0.017854 


0 .466 1 3 Q 


0.500 


0.913420 


C. 540956 


0 . 01 903 7 


0.475278 
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